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S '■ ABSTRACT 

^ ■ We investigate the dynamical evolution of compact young star clusters (CYCs) 

near the Galactic center (GC) using Aarseth's Nbody6 codes. The relatively small 
number of stars in the cluster (5,000-20,000) makes real-number N-body simulations 
for these clusters feasible on current workstations. Using Fokker-Planck (F-P) models, 
I Kim, Morris, &; Lee (1999) have made a survey of cluster lifetimes for various initial 

^—i ' conditions, and have found that clusters with a mass ^ 2 x 10^ M© evaporate in 

■ ~ 10 Myr. These results were, however, to be confirmed by N-body simulations 

QQ ■ because some extreme cluster conditions, such as strong tidal forces and a large stellar 

. mass range participating in the dynamical evolution, might violate assumptions made 

(•~^ , in F-P models. Here we find that, in most cases, the CYC lifetimes of previous F-P 

fH I calculations are 5-30 % shorter than those from the present N-body simulations. The 

comparison of projected number density profiles and stellar mass functions between 

O ' N-body simulations and HST/NICMOS observations by Figer et al. (1999) suggests 
' 

+-J ' that the current tidal radius of the Arches cluster is ~ 1.0 pc, and the following 

■ parameters for the initial conditions of that cluster: total mass of 2 X lO'' MfT, and 
J> . mass function slope for intermediate-to-massive stars of 1.75 (the Salpeter function 

^ , has 2.35). We also find that the lower stellar mass limit, the presence of primordial 

^ I binaries, the amount of initial mass segregation, and the choice of initial density profile 
(King or Plummer models) do not significantly affect the dynamical evolution of CYCs. 
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1. 



INTRODUCTION 



The Arches (GO. 121+0.017) and Quintuplet (AFGL 2004) clusters are two extraordinary 
star clusters near the Galactic center. They are very young (< 5 Myr), compact 1 pc), and 
only 20-30 parsecs away from the Galactic center (GC) in projection, while they appear to be as 
massive as the smallest Galactic globular clusters (~ 10^ Mq; Figer et al. 1999). These compact 
young clusters (CYCs) have several interesting dynamical characteristics that distinguish them 
from globular clusters: 1) CYCs have very short dynamical and half-mass relaxation timescales 
{tdyn ~ 10^~® yr and tj-h ~ 10^"''' yr, respectively); 2) CYCs are situated in strong tidal fields (the 
tidal radius of a W^Mq cluster located 30 pc from the Galactic center is ~ 1 pc); and 3) mass 
segregation may occur on a timescale shorter than the lifetimes of the most massive stars, such 
that those massive stars play an important role in the dynamical evolution of the cluster (for more 
details on timescales and the effects of tidal fields, see Kim, Morris, & Lee, 1999; hereafter KML). 

The fact that we currently observe only two such young clusters near the GC (in addition to 
the central cluster right at the GC) raises a natural question about the lifetimes of CYCs. Using 
anisotropic Fokker-Planck (F-P) models, KML surveyed lifetimes of CYCs for various initial mass 
functions (IMFs), cluster masses (M), and Galactocentric radii (Rg), and found that clusters with 



times {tev', defined here as the time by which M, the cluster mass inside the tidal radius, becomes 
5 % of its initial value) of CYCs are first due to short tdyn and t^h, and strong tidal forces, but 
the mass loss accompanying the evolution of massive stars is also responsible for shortening tev of 
clusters that last longer than ~ 3 Myr. 

F-P models are statistical models involving distribution functions. It is the statistical 
stability and fast computing time of the F-P models that the survey-type study by KML required. 
Takahashi & Portegies Zwart (1998) found good agreement between anisotropic F-P models and 
N-body simulations for globular clusters by adopting an "apocenter criterion" and an appropriate 
coefficient for the speed of star removal beyond the tidal radius (KML adopted these as well). 
However, some extreme conditions of CYCs may be inconsistent with the assumptions inherent 
to F-P models. As discussed in KML, the conditions required by F-P models, tdyn trh and 
tdyn ^ tse (tse is the Stellar evolution timescale) may be violated, especially in the core at certain 
epochs. Moreover, the active participation of a large mass range of stars in the dynamics, which is 
another peculiarity of CYCs, is difficult to realize in F-P models that embody a mass spectrum 
with a restricted number (usually 10-20) of discrete mass components. A greater number of 
components would better express the mass spectrum, but then the computation time would 



^Portegies Zwart et al. (1999) have studied with N-body models the evolution of R136, a compact, young star 
cluster in the 30 Doradus region of the Large Magellanic Cloud. This cluster resembles the CYCs near the GC in 
many ways, except for the tidal forces, which are speculated to be much weaker than those for the CYCs near the 
GC. We note that Portegies Zwart et al. (2000) are independently working on N-body simulations targeted for the 
CYCs near the GC at the moment of the submission of the present paper. 




These unparalleled, short evaporation 
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become accordingly longer. On the other hand, too small a number of components would not 
properly realize the whole mass spectrum. By fixing the number of stars in the most massive 
components, KML tried to carefully account for a relatively small number of the most massive 
stars. Yet, the results from such a treatment are to be confirmed by more realistic models, N-body 



CYCs, estimated to have ~ 10^^^ stars, are one of a few classes of systems for which 
real-number N-body simulations are feasible on current workstations. Among many virtues 
of N-body simulations, the natural realization of the mass spectrum and the tidal fields is 
particularly beneficial to the study of CYCs. These benefits allow N-body simulations to treat 
mass segregation and evaporation of stars exactly, thus providing better density profiles and mass 
spectra as a function of radius. These are photometric observables, and are partially available 
from the HST/NICMOS observations of the Arches and Quintuplet (Figer et al. 1999; the partial 
availability is owed to crowding at the cluster center and to the detection limit at the faint end of 
the luminosity function). The Arches and Quintuplet are estimated to be only ~ 2 and ~ 4 Myr 
old, respectively (Figer et al. 1999), but their currently observed structures must already have 
deviated from the initial ones due to their rapid dynamical evolution. With both real-number 
N-body simulations and HST observations of the two clusters in hand, one has a rare chance of 
deriving not only the current characteristics of these systems, but their initial conditions as well. 

In the present study, we perform a series of N-body simulations 1) to compare the lifetimes of 
CYCs with those from F-P models obtained by KML, 2) to find the initial cluster conditions that 
best match the current observations of CYCs, and 3) to test the effects of initial mass segregation 
and different initial density profiles on the dynamical evolution of the CYCs. 



For N-body simulations in the present study, we use the most recent version of Aarseth's 
N-body codes, Nbody6 (Aarseth 1999 and references therein). We modified the Nbody6 codes 
to implement the full, non-truncated tidal forces, and we adopted the prescription for stellar 
evolution used in KML for consistency (see KML for details on the stellar evolution). The upper 
stellar mass limit, m^, is 150 M©, and the lower mass limit, mi, is 0.1 M© or IMq. We set the 
cluster to initially fill the tidal radius, and remove stars outside 5 times the tidal radius. 

For initial conditions of the cluster and the Galactic tidal field, we also follow prescriptions of 
KML: We adopt single-mass King models (with a King parameter Wq = 4) and single power-law 
mass functions; The Galactic mass Mg inside a Galactocentric radius Rg is given by 



simulations. 



2. 



SIMULATION MODELS 




(1) 
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(from Genzel & Townes, 1987). Then the tidal radius of the cluster, Rt, becomes 

1/3 



/ M 



- 5.3 X 10- PC ^- J . (2) 

The equations of motions are integrated in the rotating frame, and tidal forces are calculated from 
the Galactic potential corresponding to equation (|l]). 



3. N-BODY VS. FOKKER-PLANCK 

In order to compare tev of N-body and F-P calculations, we performed N-body simulations of 
several representative models in KML that have number of stars, A^, smaller than ~ 15, 000. This 
criterion excludes nii = 0.1 Mq models with the exponent of the power-law IMF, a (defined as in 
dN oc m~'^dm), larger than or equal to 2. 

The lifetimes of 9 N-body simulations performed in the present study are shown in Table || 
along with those from F-P calculations of KML.i This set of 9 models represent different M 
(models 1 &: 5), different a (models 1, 8, 13, &: 23), different mi (models 3 &: 8; models 13 & 15), 
and different Rg (models 13, 43, & 44). We find that N-body calculations always give longer tev 
than F-P calculations, and that the fractional difference of tev between two calculations, A'tev 
(= [Nbody-FP]/Nbody), does not show any particular correlation with A^. The biggest A'tev, 
59 %, is for model 43, for which the IMF is quite flat (a = 1.75) and the cluster evaporates before 
significant stellar evolution starts. This implies that A'tev is larger for smaller a (compare to 
model 5, which has as short tev as model 43 but A'tev of only 8 %), but that the stellar evolution 
has the effect of diminishing A'tev (cf- model 44, which also has a = 1.75 but A'tev of only 16 %). 
We performed the N-body and F-P calculations for model 13 without stellar evolution, and find 
that the new tev^s are 4.5 and 2.5 Myr, respectively (model 13 is our baseline model because 
we later find that it best matches the observations; see § This indicates that the speed of 
the pure (i.e., no stellar evolution) dynamical evolution of F-P models by KML is considerably 
overestimated. 

We attribute shorter tg^'s of F-P calculations to the underestimated choice of Nu (the number 
of stars in the most massive bin; the same as when the number of mass bins, A'^bin, is 15 as 
in KML) by KML. The mass bins of F-P models by KML were chosen to make constant, so 
that the significance of the most massive bins of models with different IMFs would be the same. 



^The lifetimes from N-body simulations are subject to the statistical fluctuation which may be significant for 
models with A'' smaller than a few thousands. We find that our models have fluctuations from one run to the next 
mostly smaller than 10 %. Such relatively small fluctuations are due to our deflnition of tev (see § |l|), which avoids 
the late phase of the evolution that suffers most from statistical fluctuations. 
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The larger A^^, the smaller the characteristic mass of the most massive bin, thus the smaller the 
effective mass range (the lower and upper mass limits are not a function of N^). Since a larger 
effective mass range gives faster mass segregation and thus faster ejection of lighter stars, a smaller 
Nu results in smaller tgy. For A^bin cannot be infinitely large for a practical reason, F-P models 
have difficulties with dealing with a very large mass range, and F-P results may depend on the 
way of binning the mass with a limited number of components (bins). The choice of = 50 
in KML was arbitrary, and one may tune the F-P tev values to N-body t^y values by controlling 
this parameter. We find that a new F-P calculation for model 13 with no stellar evolution and 
Nu = 250 — 300 gives a tev value that matches the N-body results without stellar evolution, 
4.5 Myr. In conclusion, the choice of A^^^ = 50 in KML was too small, thus KML underestimated 
tev, but such underestimation was not significant for most cases because the amount of mass loss 
from stellar evolution after 2 Myr accelerates the evaporation of clusters. 

Figure || compares the volume density profiles of N-body and F-P calculations for model 13 at 
1 and 2 Myr. For this F-P calculations, = 250 and A^'bin = 10 (instead of 15 as in KML) were 
used. For Nu = 250, A'^bin larger than 10 would give negative /3 values (see KML for definition of /3; 
a negative /? gives smaller numbers of stars for lighter mass bins). At t = 1 Myr, the two density 
profiles differ only near and outside Rt (the small deviation at the core is probably due to small 
number statistics), but the F-P calculation significantly underestimates the density outside Rt- At 
t = 2 Myr, the deviation starts at smaller radii, and the power-law slope of the F-P density profile 
is steeper than the N-body slope by 0.3 to 0.4 inside Rt- Such a difference makes the cluster mass 
inside the tidal radius of the F-P calculation always somewhat smaller than that of the N-body 
calculation. Figure ^ compares the evolution of cluster mass inside the tidal radius, M, of N-body 
and F-P calculations for model 13. Stellar evolution was not included in these calculations in order 
to see the effect of tidal fields only, and A^^ = 250 and A^bin = 10 were used for the F-P calculation. 
The figure shows that the mass loss rate from the cluster differs especially during the early phase 
of evolution. We also find that the a value of the F-P calculation decreases faster with time, 
having 0.1 to 0.3 smaller (shallower mass function) values, than those of the N-body simulation 
at most times. N-body models are more realistic than F-P models, whose observables are here 
found to have non-negligible discrepancies from N-body models; therefore for comparison with the 
observed structures of the CYCs, N-body simulations should be preferred over F-P models. 

4. INITIAL CONDITIONS OF THE ARCHES 

4.1. Comparison with Observations 

In spite of the very young ages of the CYCs, the observed properties of CYCs should differ 
from their initial status because of the cluster's rapid dynamical evolution. Thus one should 
compare numerical simulations of cluster evolution with observations to infer the initial conditions 
of the cluster. As discussed in § ^, N-body simulations provide the most accurate information 
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on the dependence of the stellar distribution on stellar mass. In this section, we use our N-body 
simulations to find initial conditions of the Arches cluster from comparisons with i^ST/NICMOS 
photometric observations by Figer et al. (1999). The Quintuplet cluster has a largely dispersed 
distribution (probably because it is in the final disruption phase), so dynamical information from 
its image is very limited. 

Dynamical information from the photometric image of a cluster can be found from the surface 
density profiles and the mass distribution functions as a function of radius. When the photometry 
is complete down to the faintest stars, the overall stellar mass function, total cluster mass, and 
degree of mass segregation are obtained from the above two observables. However, the HST 
observations of the Arches were limited by crowding in the core and by the background confusion 
in the outer regions. Using the mass-K magnitude relation of Meynet et al. (1994) and adopting 
an extinction at K band of 3.1, Figer et al. (1999) estimate that their HST Arches photometry is 
complete down to m ~ 20 Mq for r > 3", or 0.12 pc, and to m ~ 8Mq for r > 5. "25, or 0.2 pc 
(the distance to the GC is assumed to be 8 kpc in the present study), and that the background 
confusion limit lies between 3 and 5Mq. For this reason, we only use stars having m > 20 Mq 
(F205W < 15 mag, where F205W is the undereddened, apparent Vega magnitude of the NICMOS 
F205W filter) for the surface number density profile, and m > SMq (F205W < 17 mag) for 
estimating the mass spectrum. The number of stars heavier than SMq outside its crowding limit, 
0.2 pc, is 232, which is not enough to give reliable information on the mass spectrum as a function 
of radius. Therefore we adopt the simplest way of measuring the mass function and the mass 
segregation: we count the number of stars in two mass bins (8 < itl/Mq < 20 and 20 < m/M©), 
each in two radius bins (0.2 < r/pc < 0.4 and 0.4 < r/pc < 0.8). In principle, one can find the 
initial conditions (or several sets of initial conditions) of a cluster by comparing these number 
counts and the density profiles of numerical simulations and observations. 

We start with IMFs that give the observed number of massive stars, which, for the Arches 
cluster, is expected to be very close to its initial value (the estimated age of the Arches, ~ 2 Myr, 
is before significant stellar evolution). The observations require the IMF to have A^(> 30 Mq) 
of 150 (30 Mq stars are complete down to r ~ 0.06 pc). The surface number density profiles for 
M > 20 Mq (S20) at 1 or 2 Myr, and the number count evolutions of N-body simulations of 
models 13, 21, and 22 are plotted in Figures ^ and |^ along with observations (the count numbers 
from observations are given in Table §). These three models have A'^(> 30 Mq) values of about 150 
but different a values (see Table ^. S20 plots in Figure |^ show that all three models agree well 
with the observation, which implies that our choice of the bulge mean density at the location of 
the cluster is appropriate because the mean bulge density determines the tidal radius of the cluster 
(this does not necessarily prove that our choices of the bulge mass profile, eq. [|| , and Rg for these 
models, 30 pc, are correct since the bulge mean density is a combination of two parameters; see 
below). The Rt values of our N-body models shown in Figure |^ are 1.0-1.2 pc. Figure ^ shows 
that the observed number counts best agree with model 13 (a = 1.75) at time t = 1-2 Myr. Model 
21 (a = 1.5) gives slightly smaller counts than observations, but considering the uncertainty in the 
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mass-magnitude relation, we cannot exclude the possibility that model 21 is applicable. On the 
other hand, model 22 (a = 2) significantly over-predicts the number of stars in lighter mass bins, 
except at t > 3 Myr, by which time the cluster would be in a disruption phase. We conclude that 
the best a value suggested by these comparisons to observations is 1.75. 

This finding may be, in fact, limited to the high mass regime only, where the comparison 
is actually made. However, for the above best-fit mass function, the stars with m > 8Mq 
already consitute more than 70 % of the total mass, thus the above a value may be considered 
representative. 

The projected distance from the GC to the Arches is ~ 20 pc, but the true distance is 
unknown. Furthermore, the bulge mass distribution at a few tens of parsecs from the GC is 
uncertain. While the infrared light distribution observed by Becklin & Neugebauer (1968) implies 
a density profile oc R^^'^ (eq. ^ exhibits this profile), the radial velocity observations in radio 
frequencies for the same region suggest a shallower density drop (oc R~^'^) between 30 and 100 pc, 
leading to a steeper enclosed-mass increase: 

M,=.8xl0'Me(^)" (3) 

(Lindqvist, Habing, &; Winnberg 1992). The tidal environment of the cluster is mostly determined 
by the mean density, Mg/Rg. We performed the best-fit model found above (a = 1.75; model 13) 
with different tidal environments: models located at Rg = 20 pc (model 19) and 50 pc (model 14) 
with our standard bulge mass distribution, equation (|l|). The latter is equivalent to the case of 
Rg = 30 pc with a bulge mass distribution of equation (|3|). We find that T,2o of models 13, 14, 
and 19 can be all nicely fit to the observed S20 at t = 1-2 Myr, and that the number counts of 
models 14 and 19 showed worse agreement to observations than model 13, but the agreement was 
still at the acceptable level. We conclude that our comparisons do not rule out the possibility of 
Rg = 20 pc and 50 pc, but favor the Rg = 30 pc case. 

Model 15 has the same initial conditions as our best-fit model, model 13, except mi = 0.1 Mq 
instead of IMq. We find that tev of model 15 is only ~ 10% longer than that of model 13 (see 
Table |^. This confirms the finding by KML that t^^ does not sensitively depend on mi when 
a < 2. This phenomenon is caused by two factors: 1) for a model with a = 1.75 and = 150 Mq, 
the mass between 0.1 Mq and IM© constitutes only 15% of the total mass, and 2) the lightest 
stars are rapidly ejected from the cluster during the early phases due to the large m^/mi ratio. 
However, the number counts of model 15 agree well with observations only after t = 2.4 Myr, by 
which time the Arches cluster would be in the disruption phase and would have a structure as 
dispersed as the Quintuplet. Thus the observations more support our model with mi = 1 M© than 

0.1 Mq. 

One of the most commonly used IMFs for the Galactic disk, the Scalo mass function (Scalo 
1986), has a = 2.7 for m > IM0, implying that the best-fit a value obtained here for the 
Arches cluster is considerably flatter than that for the disk. This fact, together with mi = 1 Mq 
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being favored over mi = 0.1 M©, appears to support the arguments by Morris (1993) that the 
non-standard star formation environment near the GC may lead to an IMF skewed toward 
relatively massive stars and having an elevated lower mass cutoff. 



4.2. Primordial Binaries 

So far, we have not considered primordial binaries. Camera 2 of the HST/NICMOS instrument 
(used by Figer et al. 1999) has an angular resolution of about 1 pixel size of the detector (0.075"), 
which is ~ 600 AU at the distance of the GC. Thus binary systems with a semi-major axis smaller 
than few hundred AU are not resolved, and the presence of a significant number of primordial 
binaries in the cluster may affect our number count analysis above by decreasing number counts 
and/or by moving a primary star to a more massive mass bin (the number of dynamical binaries, 
i.e., binaries formed through close encounters, at a given moment is only a few, if any, for the whole 
cluster lifetime). To see the effects of primordial binaries on the number counts, we performed two 
simulations of model 13 with primordial binary fractions, fbin, of 25 % and 50 %. The fraction is 
defined as ^ 

where Nt,in and Nging are the numbers of binary systems and single stars, respectively. Thus the 
percentages of stars in binary systems for fi^n = 25 % and 50 % are 40 % and 67 %, respectively. 
These relatively 'moderate' fbin values, compared to > 50 % used in N-body simulations 
for open clusters (Kroupa, Petr, & McCaughrean 1999, for example) may be justified by an 
argument by Durisen &: Sterzik (1994) that binary formation from fragmentation of collapsing 
and rotating clouds or from a gravitational instability of massive protostellar disks is more 
likely in low-temperature clouds (the central molecular zone in the inner few hundred pc of the 
Galactic bulge has significantly elevated temperatures, ~ 70 K). The initial companion-mass-ratio 
distribution is obtained by random pairing of stars, and the initial eccentricity distribution is 
assumed to be thermally relaxed. For the initial period distribution, we adopt equation (8) of 
Kroupa (1995a), which approximates the distribution of binary systems in the Galactic disk. 
This initial distribution is evolved prior to the start of the N-body integration to account for the 
"pre-main-sequence eigenevolution" (the evolution in orbital parameters due to internal processes 
such as tidal circularization during the pre-main-sequence; see Kroupa 1995a for details). 

We find that the tev values of fbin = 25 % and 50 % models are slightly (10-20 %) longer 
than that of fbin = % model. This insensitivity of 

tev Oil fbin was noted by McMillan & Hut 
(1994), and Kroupa (1995b). Kroupa et al. (1999) find that their model with fbin = 100 % evolves 
slightly more slowly than that with fbin = 60 % and interpret it due to cooling by disruption of 
wide binaries. The slightly longer tgv's of our models with primordial binaries may be interpreted 
in the same way, but the reduced number of effective point sources due to binarity may be 
another possible explanation (the rate of relaxation is proportional to Nm^, where m is the stellar 
mass; when companions are picked out of a pool with a large mass range, the binarity does not 
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significantly increase the m of an 'effective point mass', wliile it considerably decreases of the 
effective point masses). 

The number counts with a consideration of the angular resolution of the camera for 
fhin = 25 % and 50 % models are compared to those of the fhin = % model in Figure ^. Here, 
the mass of a binary system with a semi-major axis smaller than the angular resolution is obtained 
from the total luminosity of the binary system, and the mass-luminosity relation adopted in Figer 
et al. (1999) is used for this calculation. While the number counts of the fbin = 25 % model are 
only slightly smaller than the fbin = % model, those of the fbin = 50 % model are 10-40 % 
smaller. However, these relatively smaller number counts of models with primordial binaries still 
agree with observed number counts at an acceptable level. Furthermore, the difference of number 
counts between the fbin = % and 50 % models is larger for the more massive bin, implying that 
the initial a should be smaller than that of the model tried here (1.75) in order for the fbin = 50 % 
model to better match the observations. In conclusion, the presence of primordial binaries will not 
significantly change our findings above on the best-fit initial conditions of the Arches cluster, and 
our results indicate that the a value found above for the case with no primordial binaries, 1.75, is 
an upper limit. 

5. INITIAL CLUSTER STRUCTURE 

5.1. Initial Mass Segregation 

The prior sections of the present study and KML assumed, for simplicity, the same initial 
density profiles for all stellar masses, i.e., no initial mass segregation. This assumption is, however, 
not based on any observational or theoretical evidence. Interactions between stars during the 
star formation process determine the degree of initial mass segregation, which is therefore an 
important piece of information for the theory of cluster formation. However, different models lead 
to totally different predictions for initial segregation. The model by Podsiadlowski &: Price (1992), 
in which favorable stellar masses are determined by the ratio of the timescales for protostellar 
collisions and of gas infall onto the protostars, predicts more massive stars at larger radii. On the 
other hand, segregation of more massive stars in the core is predicted in some models such as the 
one by Murray & Lin (1996), where encounters between cloudlets increase the protostellar masses, 
and the one by Bonnell et al. (1997), where the deeper potential in the core causes stars there to 
accrete more circumstellar material. Here we attempt both cases, one with heavy stars initially 
more prevalent at the core, the other with heavy stars more in the envelope. 

The initial density and velocity profiles of model 13 are modified to include the initial mass 
segregation: model 31 initially has a King profile with Wq = 2 for the lightest stars and Wq = 6 
for the heaviest stars (more heavy stars in the envelope) while model 32 has the opposite Wq 
values (more heavy stars at the core). The intermediate-mass stars have interpolated Wq values 
depending on the logarithm of their masses. We find that models 13, 31, & 32 show very similar 
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profiles at 2 Myr. In spite of different initial conditions, the number counts of these three 
models also exhibit similar evolution after t = 1 Myr (see Fig. |6|). It appears that the relaxation 
processes are rapid enough to erase the memory of the initial mass segregation in less than 1 Myr. 
Thus unless the initial segregation is more severe than models tested here, the best-fit initial 
conditions found in § ^ are robust against the initial segregation. 

5.2. King vs. Plummer Models 

KML showed that the global evolution (such as M & i?^) of clusters initially having King 
profiles of Wq = 1-7 does not depend on Wq values. This is because King profiles with different 
Wq values mainly differ at the core while the global evolution of the cluster is mostly determined 
by the properties in the cluster envelope. Here we make two N-body simulations initially having 
Plummer models to see the effects of initial conditions different from King models. Plummer 
models have a density profile of 

POC (i?2 + i?2)-5/2^ 

which has no density drop analogous to tidal cutoff in King models (i?c is the core radius). Thus 
the density profiles of Plummer models have shallower decrease than King models near Rt, and 
have to be artificially cut at Rf. The only parameter that determines the profile is the ratio of 
Rt to Rc, Cp. Here we perform two simulations analogous to model 13 with Cp = 3 and Cp = 6. 
The Plummer model with Cp = 6 has core density and half-mass radius comparable to those of the 
King model with Wq = 4. We find that the Plummer models give tev values only slightly longer 
than that of the corresponding King model (by less than 20 %), and S20 and the number counts 
of the Plummer models are very similar to those of the King model. Therefore, the findings in § ^ 
may also apply to models with Plummer initial conditions with Cp = 3-6. 

6. SUMMARY 

Using Aarseth's Nbody6 codes, we have studied the dynamical evolution of CYCs near the 
GC. First, we confirm the results of KML that clusters with a mass ^ 2 x 10^ Mq evaporate in 
~ 10 Myr, but find that the F-P calculations by KML underestimated the lifetimes of CYCs by 
5 to 30 % in most cases. This discrepancy is due to the adoption by KML of values of Nu which 
are too small, and we find that Nu of 250-300 would be appropriate although such a large number 
would bring the effctive mass of the largest mass bin down to a value at which we would lose some 
of the effects associated with the presence of the most massive stars. Without stellar evolution, 
the above t^v discrepancy would be more considerable, i.e. the mass loss from stellar evolution 
starting slightly after t = 2 Myr significantly accelerates the evaporation of the clusters. The 
volume density profiles from F-P calculations are steeper than those from N-body simulations, 
and especially, the F-P densities outside Rt are significantly underestimated. 
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By comparing the surface number density profiles and number counts in 2 mass bins at 2 
radius bins between N-body simulations and il/^ST/NICMOS photometry of the Arches cluster, we 
find the following best-fit initial conditions: M = 2 x 10"^ Mq and a = 1.75. The relation between 
M and a is constrained by the observed number of stars heavier than 30 M©, 150. Larger or 
smaller a values than 1.75 give less satisfactory number count fits to observations. The presence 
of primordial binaries favors slightly smaller a values than 1.75. The fit of S20 to observations 
indicates i?^ ~ 1 pc, but Rg is not well constrained. The mean density of the bulge inside the 
Arches cluster is suggested to be 5-30 x IO^Mqpc"^. Also, we confirm the finding of KML that 
for clusters with a < 2, t^v of a cluster with m/ = 0.1 is very close to that of a cluster with 
mi = 1M0. However, the number count plots seem to more support m; = IMq. These findings 
are in agreement with the arguments by Morris (1993) that the non-standard star formation 
environment near the GC may lead to an IMF skewed toward relatively massive stars and having 
an elevated lower mass cutoff. 

The global evolution and tgv of CYCs are nearly independent of the initial mass segregation 
and the choice of initial density profile. Clusters with initial mass segregation and reverse 
segregation (more heavy stars in the outer region) show very similar evolution except in the very 
early phase (< 0.5 Myr). Plummer initial models with Cp = 3 and 6 exhibit indistinguishable 
evolution from a King initial model. Therefore, our findings regarding the initial conditions of the 
Arches cluster appear to be insensitive to the initial mass distribution. 

As discussed in KML, the effect of the gas left over from cluster formation may be important 
to the early dynamical evolution of CYCs because the remnant gas is thought to be blown away 
from the cluster in the early phases by strong stellar winds or supernova explosions, and such 
abrupt disappearance may significantly change the potential of the cluster in a short amount of 
time. Consequently, the role of the remnant gas in early cluster evolution will be the subject of a 
future study using N-body simulations. 

S.S.K. is deeply grateful to Sverre Aarseth for generously providing us with his Nbody6 codes 
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Table 1. Comparison of tgv between N-body & F-P Calculations 







M 




mi 




Rg 




(Myr) 


A't 


Model 


KML 


(Mo) 


a 


(Mo) 


N 


(pc) 


N-body F-P 


(%) 


1 


142 


2 X 10"^ 


2.35 


1.0 


6270 


30 


5.5 


4.0 


27 


3 


115 


2 X 10'' 


1.50 


0.1 


5164 


30 


2.8 


2.2 


21 


5 


141 


5 X 10^ 


2.35 


1.0 


1567 


30 


2.4 


2.2 


8 


8 


113 


2 X 10^ 


1.50 


1.0 


1633 


30 


2.7 


2.3 


15 


13 




2 X 10^ 


1.75 


1.0 


2605 


30 


2.8 


2.4 


14 


15 




2 X 10^ 


1.75 


0.1 


12706 


30 


3.0 


2.6 


13 


23 


101 


2 X 10^ 


2.00 


1.0 


6270 


30 


3.6 


2.7 


25 


43 




2 X lO'' 


1.75 


1.0 


2605 


10 


2.2 


0.9 


59 


44 




2 X 10^ 


1.75 


1.0 


2605 


100 


4.5 


3.8 


16 



Note. — F-P values are from KML except for models with a = 1.75, which 
were additionally calculated for the present stiidy. The coefficient for the speed of 
star removal, a^sc, of 2 is adopted for all F-P simulations here. Column KML is 
for the model numbers from the paper KML. A'tev is the fractional difference, |FP- 
Nbodyl/Nbody. 



Table 2. Number Counts from the Observation of the Arches Cluster 



Mass Bin 


Radius Bin 


8 < m/M© 


< 20 20 < 


ui/Mq 


0.4 < r/pc < 0.8 


64 




46 


0.2 < r/pc < 0.4 


66 




56 
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Table 3. N-body Simulations for the Arches' Initial Conditions 





M 




mi 


Rg 




Model 


(Mo) 


a 


(Mo) 


(pc) 


(Myr) 


13 


2.0 X 10'' 


1.75 


1.0 


30 


2.8 


14 


2.0 X 10^ 


1.75 


1.0 


50 


4.1 


15 


2.0 X 10^ 


1.75 


0.1 


30 


3.0 


19 


2.0 X 10^ 


1.75 


1.0 


20 


2.4 


21 


1.6 X lO'' 


1.50 


1.0 


30 


2.6 


22 


2.8 X 10^ 


2.00 


1.0 


30 


3.7 
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Fig. 1. — Volume density profiles of N-body (thin lines) and F-P (thick lines) calculations for model 
13 at t = 1 Myr {upper panel) and 2 Myr (lower panel). The locations of tidal radii are marked 
with arrows. Nu = 250 and AT^iu = 10 were used for the F-P calculation. 
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Fig. 2. — Evolution of cluster mass inside the tidal radius, M, of N-body (thin lines) and F-P (thick 
lines) calculations for model 13. M's are normalized to their initial values. Stellar evolution is not 
included in these calculations, and A^^ = 250 and iVbin = 10 were used for the F-P calculation. 
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Fig. 3. — Surface number density profiles of stars heavier than 20 Mq, E20, from our N-body 
models 21 {a = 1.5; triangles), 13 (a = 1.75; squares), & 22 (a = 2.0; pentagons), and NICMOS 
observations of the Arches by Figer et al. (1999; solid circles connected with a line). N-body results 
shown are the ones at t = 2 Myr for models 13 &; 22, and t = 1 Myr for model 21 (model 21 shows 
a worse fit for the core region at t = 2 Myr). The tidal radii of N-body models 21, 13, &; 22 at 
the epoch shown here are 1.0 pc, 1.0 pc, & 1.2 pc, respectively. The observations are limited by 
crowding at the core and by background confusion near and outside the tidal radius. 1-a Poisson 
errors for observations are indicated by vertical bars. 
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Fig. 4. — Evolution of number counts of stars with 8 < m/M© < 20 (thin lines) and stars with 
20 < ui/Mq (thick lines) in the 0.2 < r/pc < 0.4 region {lower panels) and the 0.4 < r/pc < 0.8 
region (upper panels), from our N-body models 21 (a = 1.5; solid lines in left panels), 13 (a = 1.75; 
solid lines in middle panels), &l 22 (a = 2.0; solid lines in right panels), and NICMOS observations 
of the Arches by Figer et al. (1999; dashed lines). The dashed lines cover the estimated age of the 
Arches cluster, 2 it 1 Myr. 
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Fig. 5. — Evolution of number counts of stars with 8 < m/M© < 20 {thin lines) and stars with 
20 < m/MQ {thick lines) in the 0.2 < r/pc < 0.4 region {lower panels) and the 0.4 < r/pc < 0.8 

region {upper panels), from our N-body model 13 with = 25 % {dashed lines in left panels) 
and model 13 with f^in = 50 % {dashed lines in right panels). The angular resolution of the 
HST/N1CM.0S camera 2, ~ 0.075", is considered in counting apparent numbers and estimating 
apparent masses (see text for details). Model 13 with /f,j„ = % {solid lines in both panels) is also 
plotted for comparison. 
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Fig. 6. — Evolution of number counts of stars with 8 < m/M© < 20 {thin lines) & stars with 
20 < m/M(7) (thick lines) in the 0.2 < r/pc < 0.4 region (lower panel) and the 0.4 < r/pc < 0.8 
region (upper panel), from our N-body models 13 (no initial segregation; solid lines), 31 (initially 
more heavy stars in the evelope; dashed lines), k. 32 (initially more heavy stars in the core; dotted 
lines) . 



